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Abstract 

We report on a significant failure of the local density approximation (LDA) and the generalized 
gradient approximation (GGA) to reproduce the phase stability and thermodynamics of mixed- 
valence Li x FeP04 compounds. Experimentally, Li x FeP04 compositions (0 < x < 1) are known 
to be unstable and phase separate into LiFeP04 and FeP04. However, first-principles calcula- 
tions with LDA/GGA yield energetically favorable intermediate compounds and hence no phase 
separation. This qualitative failure of LDA/GGA seems to have its origin in the LDA/GGA self- 
interaction which delocalizes charge over the mixed-valence Fe ions, and is corrected by explicitly 
considering correlation effects in this material. This is demonstrated with LDA+U calculations 
which correctly predict phase separation in Li x FeP04 for U — J > 3.5eV. The origin of the desta- 
bilization of intermediate compounds is identified as electron localization and charge ordering at 
different iron sites. Introduction of correlation also yields more accurate electrochemical reaction 
energies between FeP04/Li x FeP04 and Li/Li + electrodes. 

PACS numbers: 71.15.Mb, 71.27. +a, 91.60.Ed 
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First-principles calculations employing density functional theory (DFT) have proven to 
be a powerful method in understanding the thermodynamic, structural and electronic prop- 
erties of a large class of materials. The density functional is not known exactly, and is 
usually modeled within the Local Density Approximation (LDA) or Generalized Gradient 
Approximation (GGA). For many systems LDA or GGA gives remarkably good agreement 
with experiments, which has made these techniques valuable tools to predict the behavior of 
materials [lj]. However, the self-interaction in LDA/GGA tends to delocalize electrons too 
much, and as such these methods are unable to capture precisely the Coulomb correlation 
effects in correlated electron systems like transition metal oxides. The resultingfailure to 
predict many transition metal oxides as insulators has been well documented |2j. In this 
paper we show by means of olivine-type LL J FeP04 that the tendency for LDA/GGA to de- 
localize the ci-electrons in mixed-valence transition metal oxides also leads to a qualitative 
failure in predicting miscibility and phase stability by a surprisingly large magnitude. The 
role Coulombic correlations play in phase stability will be qualitatively probed. 

LiFeP04, a naturally occurring mineral, has attracted much attention recently, as its 
superb thermal safety, non-toxicity and low cost make it the most likely candidate for 
rechargeable Li-batteries electrodes in large applications such as electric and hybrid ve- 
Mel. JflBflfla. Ina b a ttery , Uthinn, is elect—call, and r eve relbly cycled in 
and out of the LiFeP04 material. As a result, the pseudo-binary FeP04 - LiFeP04 phase 
diagram, critical for the material's behavior as an electrode, has been well characterized 
exp er iment ally. 

Olivine-type LiFeP04 and the de-lithiated structure FeP04 have an orthorhombic unit 
cell with four Formula Units (FU) and space group Pnma (see Fig. Q). The olivine-type 
structure contains a distorted hexagonal close-packing of oxygen anions, with three types of 
cations occupying the interstitial sites: 1) corner-sharing FeOg octahedra which are nearly 
coplanar to form a distorted 2-d square lattice perpendicular to the a axis, 2) edge-sharing 
LiOe octahedra aligned in parallel chains along the b axis, and 3) tetrahedral PO4 groups 
connecting neighboring planes or arrays. Electrochemical experiments and X-ray diffraction 
measurements have confirmed that no intermediate compound Lf r FeP04 exists between 
FeP0 4 and LiFeP0 4 0, Q, so that its phase diagram consists of a wide two-phase region 
with limited solubility on both the FeP04 and LiFeP04 sides. The ma gne tic structure 
of LiFeP04 and FeP04 was determined from neutron diffraction data fcj |.10l. Below the 
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Neel temperature Tn = 50K [9 J and 125K [10], respectively, the iron spins align in an 
antiferromagnetic (AFM) array, induced by Fe-O-Fe superexchange interactions between 
neighboring iron atoms. 

The objective of this paper is to investigate the stability of compounds between the 
composition FePC>4 and LiFePC>4 and demonstrate that Coulomb correlations are essential 
in reproducing the absence of intermediate compounds. Different Li arrangements with 4 
formula units are considered in the primitive cell. All possible symmetry- distinct decorations 
of the 4 Li sites gives seven structures, including two end members (x=0, 1), one structure 
at each of x=0.25 and 0.75, and 3 at x=0.5, here named 0.5a, 0.5b and 0.5c. The structures 
0.5a, 0.5b and 0.5c have Li remaining at sites 1 and 2, 1 and 3, and 1 and 4, respectively (see 
Table H} . All the five intermediate structures have lower symmetry than the end members, 
and are monoclinic or triclinic. Total energy calculations were performed for the seven 
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TABLE I: Fractional positions of the four Li and four Fe atoms within the unit cell. 

structures in GGA (or LDA where explicitly stated) with the projector-augmented-wave 
(PAW) method [ill . Il^ as implemented in the Vienna Ab-initio Simulation Package (l^ . 
An energy cut-off of 500 eV and appropriate k-point mesh were chosen so that the total 
ground state energy is converged within 3meV per FU. All the atoms and cell parameters 
are fully relaxed at each structure. For x=0.25 and 0.75 the remaining S2 point group 
symmetry has to be removed by imposing different initial magnetization on the irons to get 
the electronic ground state (see below). The results in this paper represent the ferromagnetic 
(FM) spin-polarized configurations unless stated explicitly. Although the magnetic ground 
state of LiFePC-4 and FeP0 4 is AFM Q, Q, the difference in FM and AFM formation 
energies (defined below) is a few meV/FU in most cases, not exceeding 12 meV, and does 
not affect the qualitative analysis, which is clearer in the FM configuration. Iron is found 
to be always in the high-spin state, with the five majority spin 3<i-orbitals occupied. 
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Here we define AE(x), the formation energy per FU of Li r FeP04 as 

AE(x) = E(x) -(xE(x = l) + (l-x)E(x = 0)) (1) 

where E(x) is the ground state total energy per FU for the structure with lithium concentra- 
tion x. A negative formation energy means compound formation is energetically favorable. 
In order for phase separation to occur at room temperature, all intermediate structures 
should have positive formation energy, large enough to overcome the potential entropy gain 
in mixing. LDA results of AE(x) for all five structures are negative. Although GGA slightly 
increases the formation energy, the prediction remains qualitatively in disagreement with ex- 
periment. 
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TABLE II: LDA and GGA formation energy (meV/FU) at different Li concentrations. 

Given that the true formation energies should all be positive, these errors are large and 
somewhat surprising, since formation energies are properly weighted energy differences be- 
tween similar structures, and as such usually benefit from significant error cancellations. For 
example, in many binary alloys formation energies are only 100~200meV/atom in magni- 
tude, and hence large errors such as those found here would make them completely unreliable, 
which, based on the good agreement of many LDA/GGA studies with experiment, is not 
the case jl^j ]. 

To investigate whether Coulombic on-site effects could be related to this substantial fail- 
ure of LDA/GGA we carried out rotationally invariant LDA+U (GGA+U, more accurately) 
jl^l ] calculations. The essence of the method can be summarized by the expression for the 
total energy 

E hDA+u [p, n] = E hDA [p\ + E Hnh [h] - E dc [n] = E hDA [p] + E v [h] (2) 

where p denotes the charge density and n is the iron on-site 3d occupation matrix. The 
Hatree-Fock like interaction -Efiub from the Hubbard model replaces the double counting (dc) 
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term E^c representing the LDA on-site interaction. The U correction term Eu = -EhiA — -Eke 
is defined by Eq. |21 However E^ c is not uniquely defined, and here we consider three common 

n n 

approaches [lq . The "around mean field" dc functional [17j (del) yields low-spin iron, in 
disagreement with experiment 10]. This is not surprising since del usually works poorly in 
strongly correlated systems. We then compared formation energies with the dc functional 
defined in (dc2) and with its spherically averaged version [18] (dc3). The latter reads 

Edc(n) = — - — Trn(Tm — 1) = —Trh(Tih — 1), (3) 

Eu(n) = ^^Tr (n(l - n)) = y Tr (n(l - n)) , (4) 

where we have defined U' = U — J. The formation energy with dc2 is very insensitive to a 
large range of J (0- 2eV) when U' is fixed, and agrees with dc3 results within 10 meV for 
U' > 2eV. Therefore, we will use dc3, where there is only one effective parameter, V . We 
evaluate all results as function of U', spanning the range from to 5.5eV. When calculating 
formation energies for a given V, we assume U' to be the same for all structures. The choice 
of U' is a source of uncertainty in LDA+[7 calculations. However, we present the results 
as a function of U' and will argue that the correct physics is obtained within a reasonable 
range of V. 

In Fig. El formation energies at different U' are shown as a function of Li concentration x. 
At each concentration AE increases with U' and becomes positive at intermediate U'(~2.5- 
3.5eV). The formation energies saturate to a nearly constant value around U' ~3.5-4.5eV. 
The effect of the Ejj term is to drive the Fe-3d orbital occupation numbers to integer (0 
or 1) values. As a result, the Fe ions tend to have integral occupancy even in the partially 
lithiated structures, and charge ordering occurs: we see distinct Fe 3+ and Fe 2+ in DFT+U 
instead of the uniform Fe^ 3_:c ^ + seen in LDA/GGA. For low U' values (U' < leV) the four 
Fe ions in the unit cell have similar 3d electron occupancy and Fe-0 bond lengths for all the 
intermediate structures. Therefore, little charge ordering occurs in this limit, even though 
the Fe ions occupy symmetrically distinct positions. We will call these Fe cations (3-x) + 
like. They are stable with respect to small perturbations in initial charge distribution. In 
the high limit of C/'(>3.5 or 4.5 eV) there are 2 types of Fe ions, one very similar to those 
in FePC>4 (which we call Fe 3+ like) and the other similar to those in LiFeP04 (called Fe 2+ 
like). The designation 3(2)+ is only meaningful in that the Fe ions are similar to those 
in FeP04(LiFeP04). The Fe-0 hybridization gives them less than nominal charge. For 



x=0.25(0.75) calculations imposing the symmetry of the structure on the charge density 
leads to two 3(2)+ like and two 2.5+ like Fe ions. Only when symmetry is broken does a 
lower energy state with three 3(2)+ like and one 2(3)+ like ions form. In these structures the 
charge density has lower symmetry than what would be expected from the ionic positions 
and, hence, charge ordering occurs. As the analysis for all five structures is similar we choose 
x=0.5 typical intermediate structure for furthur discussion. 

In Fig.|3]A.E(x=0.5a) is shown as a function of U'. We investigated AFM spin configura- 
tions in x=0, 0.5a and 1 and found them to give only slightly lower total energies. The AFM 
AE (dotted line) is almost equivalent to the FM one with charge ordering (solid line). We 
also studied a 'restricted' FM system at x=0.5a where all four Fe ions have the same initial 
magnetization, ending up 2.5+ like. Charge ordering is absent in this metastable state, 
which has higher total energy than the charge-ordered ground state. From Fig. |3] we can 
compare AE with and without charge ordering. Note that the curve with charge ordering 
levels off for U > 4.5eV, which is explained below. 

To study quantitatively the change in formation energies and electron distribution as U' 
is increased, the contributions to AE are separated into the LDA energy, A^lda, and the 
correction term, AEjj, with definitions analogous to AE in Eq. ^ The occupancy of the 
most occupied of the five minority-spin 3c?-orbitals of iron is displayed in the lower part of 
Fig. HJ This orbital is most relevant because its occupation makes the difference between 
Fe 3+ and Fe 2+ . When charge ordering is absent, the occupation number does not change 
much with U' and stays near 0.5, as expected of a 2.5+ like Fe cation. In contrast, the curves 
in the charge-ordered state separate beyond U' ~ leV, with half of the ions becoming 2+ 
like and the other half 3+ like. These occupancies can help to explain AE in the upper part 
of the diagram. When charge ordering is absent (dotted lines) the four Fe cations in the 
x=0.5a structure are equally affected by U' in terms of 3d occupation, as they are in x=0 and 
1, and the changes in Tr(n(l — h)) in Li .5 a FePO 4 are canceled by the weighted average of 
those in x=0 and 1 structures. As a result, the correction term AEu is almost proportional 
to V, explaining its linear behavior in Fig. |3J When the symmetry is sufficiently broken, 
AEu w ih make Fe-3<i charge density order so as to create, as much as possible, orbitals 
with integer occupation. This comes at a cost to AEldAj which changes from large negative 
values at U' = to positive values. We see two possible reasons why A£lda increases when 
charge ordering occurs. Localization of the minority spin electrons into half of the Fe sites 
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as Fe 2+ obviously leads to an increase in kinetic energy. Additionally, since Fe 2+ and Fe 3+ 
have different Fe-0 bond lengths, their coexistence in one structure comes with a penalty 
in elastic energy. As the increase, relative to the LDA values, in AEjj is much smaller than 
in Ai?L D A, the latter can be identified as the cause of phase separation. 

The ground state electronic structure is also affected. The x=0.5a compound is insulating 
when charge ordering occurs in LDA+[7, while it is metallic in LDA or LDA+U without 
charge ordering. The end members are insulating in both LDA and LDA+£7. 

A weakness of the LDA+[7 method is that U is an external parameter, and some justi- 
fication for the choice of it is required. Considering a realistic J = leV we find phase 
separation in the LLTePC^ system for U > 3.5eV+ J=4.5eV. Above this cutoff the forma- 
tion energies and orbital occupancies become less sensitive to U. The value of U for these 
systems is likely to be even higher than this cutoff. A recent ab-initio computation of U — J 
in the related Fe2SiC>4 fayalite system suggests a value of 4.5eV for iron [19(. 

Another way to determine a physical value of U is to compare the calculated and exper- 
imental reaction energy of FePC>4 and Li to form LiFePC>4. 

-^reaction = (-^Li + -EW04 — -^LiFePO-J- (5) 

In this reaction Li + is inserted into the FePC>4 host and an electron is added to the in- 
states, reducing Fe 3+ to Fe 2+ . Since the electron addition energy for Fe 3+ is a significant 
component of this reaction energy, the result will depend on the value of U. Experimentally, 
this energy can be measured very accurately, as it is the equilibrium electrical potential 
between Li^FePC^ and Li-metal electrodes in a Li-electrolyte. In Fig. |5] the calculated 
potential is plotted as a function of U with FM and AFM spin configurations, respectively. 
The experimental voltage of 3.5V Q is reached at U — J ~ 4.2eV. 

We have further confirmed that the positive formation energies obtained in Fig. |3] are 
not an artifact of using a single unit cell by calculating the energy of four other structures 
(x=0.25 or 0.75) with a doubled unit cell. We found all these formation energies to be within 
±10meV of the results shown in Fig. 01 Positive formation energies in GGA was recently 
confirmed in Ref. j^J. 

In summary, we find that both LDA and GGA qualitatively fail to reproduce the ex- 
perimentally observed phase stability and mixing energetics in the Li a .FeP04 system. For 
U — J > 3.5eV, LDA+[7 calculations give positive AE, in agreement with experiments. 
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Hence, we speculate that the experimentally observed phase separation is due to the cost 
in kinetic and elastic energies when Fe 2+ and Fe 3+ coexist in Li x FeP04 structures. This 
physics is not well captured by LDA/GGA, as the self-interaction causes a derealization of 
the (i-electrons, resulting in electronically identical Fe ions. As a result, there is no phase 
separation in LDA/GGA, in clear disagreement with experiment. 
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Dr. E. Wu for his help in computation. 



* Electronic address: fzhou@mit.edu 
t Electronic address: gceder@mit.edu 



[1] For review, see G. Ceder, Science 280, 1099 (1998). 

[2] For reviews, see Strong Coulomb Correlations in Electronic Structure calculations, edited by 
V. I. Anisimov, Advances in Condensed Matter Science (Gordon and Breach, 2000). 

[3] A. K. Padhi, K. S. Nanjundaswamy, and J. B. Goodenough, J. Electrochem. Soc, 144, 1188 
(1997). 

[4] A. K. Padhi, K. S. Nanjundaswamy, C. Masquelier, S. Okada, and J. B. Goodenough, J. 

Electrochem. Soc, 144, 1609 (1997). 
[5] A. S. Andersson, J. O. Thomas, B. Kalska and L. Haggstrom, Electrochem. Solid State Lett, 

3, 66 (2000). 

[6] A. Yamada, S. C, Chung, and K. Hinokuma, J. Electrochem. Soc, 148, A224 (2001) 
[7] H. Huang., S-C. Yin, and L. F. Nazar, Electrochem. Solid State Lett. 4, A170 (2001). 
[8] P.P. Prosini, D. Zane and M. Pasquali, Electrochem. Acta 46, 3517 (2001). 
[9] R. P. Santoro and R. E. Newnham, Acta Cryst. 22, 344 (1967). 
[10] G. Rousse, J. Rodriguez-Carvajal, S. Patoux, and C. Masquelier, Chem. Mater. 15, 4082 
(2003). 

[11] PE. Blochl, Phys. Rev. B 50, 17953 (1994). 

[12] G. Kresse, and D. Joubert, Phys. Rev. B 59, 1758 (1999). 

[13] G. Kresse and J. Furthmuller, Phys. Rev. B 54, 11169 (1996); Comput. Mater. Sci. 6, 15 
(1996). 



8 



FIG. 1: Structure of LiFeP04with cation polyhedra. 



[14] for review see D. de Fontaine, Solid State Physics, edited by H. Ehrenreich and D. Turnbull 

(Academic Press, 1994). 
[15] A. I. Liechtenstein, V. I. Anisimov, J. Zaanen, Phys. Rev. B 52, R5467 (1995). 
[16] see A. G. Petukhov, I. I. Mazin, L. Chioncel and A. I. Lichtenstein, Phys. Rev. B 67, 153106 

(2003) and references therein. 
[17] V. I. Anisimov et al, Phys. Rev. B 48, 16929 (1993). 

[18] S. L. Dudarev, G. A. Botton, S. Y. Savrasov, C. J. Humphreys and A. P. Sutton, Phys. Rev. 
B 57, 1505 (1998). 

[19] M. Cococcioni, A. Dal Corso, and S. de Gironcoli, Phys. Rev. B 67, 094106 (2003). 
[20] P. Wang, and N. A. W. Holzwarth, Phys. Rev. B 68, 165107 (2003). 



9 




FIG. 2: Formation energy of Li^FePO^ at different x and U' values. Points at x= 0.5 correspond 
to struc 
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FIG. 4: Upper part: LDA (triangle) and U correction term (diamond) contributions to Ai?(x=0.5a) 
vs. U'. Solid/dotted lines indicate presence/absence of charge ordering. Lower part: occupancy of 
the most occupied minority-spin orbital versus U', for Fe 2+ (solid line) and 3+ (dashed line) in 
the charge-ordered state and for 2.5+ (dotted line) in the state without charge ordering. 
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FIG. 5: Reaction energy in Eq. [5]per FU versus U' with FM and AFM configurations, respectively. 
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